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STEADY  RIGID-BODY  MOTIONS 
IN  A  CENTRAL  GRAVITATIONAL  FIELD 

Li-Sheng  Wang  *  J.H.  Maddocks*  P.S.  Krishnaprasad  f 


ABSTRACT. 

In  recent  work,  the  exact  dynamic  equations  for  the  motion  of  a  finite  rigid 
body  in  a  central  gravitational  field  were  shown  to  be  of  Hamiltonian  form 
with  a  noncanonical  structure.  In  this  paper,  the  notion  of  relative  equi¬ 
librium  is  introduced  based  upon  this  exact  model.  In  relative  equilibrium, 
the  orbit  of  the  center  of  mass  of  the  rigid  body  is  a  circle,  but  the  center 
of  attraction  may  or  may  not  lie  at  the  center  of  the  orbit.  This  feature 
is  used  to  classify  great-circle  and  non- great- circle  orbits.  The  existence  of 
non-great-circle  relative  equilibria  for  the  exact  model  is  proved  from  var¬ 
ious  variational  principles.  While  the  orbital  offset  of  the  non-great-circle 
solutions  is  necessarily  small,  a  numerical  study  reveals  that  there  can  be 
significant  changes  in  orientation  away  from  the  classic  Lagrange  relative 
equilibria,  which  are  solutions  of  an  approximate  model. 


1.  Introduction 

This  paper  describes  some  results  concerning  the  relative  equilibria,  or  steady  circular 
orbits,  of  a  rigid  body  of  finite  extent  moving  under  the  effects  of  an  inverse  square 
gravitational  force  field.  In  different  parameter  regimes,  the  analysis  applies  equally  well 
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to  an  artificial  satellite  orbiting  the  Earth,  to  a  moon  orbiting  a  planet,  or  to  a  planet 
orbiting  a  sun.  In  any  of  these  cases  the  approximate  problem  at  hand  is  that  of  the 
motion  of  a  two-body  system.  We  shall  further  assume  that  a  restricted  two  body  problem 
is  appropriate,  i.e.  one  of  the  bodies  is  extremely  massive  in  comparison  to  the  other,  so 
that  the  motion  of  the  heavier  body  is  unaffected  by  the  motion  of  the  lighter  one. 

Of  course  if  the  two  bodies  are  further  approximated  as  two  point  masses  the 
restricted  two  body  problem  is  completely  integrable  and  is  completely  solved.  The  points 
of  interest  that  will  be  addressed  here  arise  due  to  effects  arising  from  the  finite  extent  of  the 
bodies.  This  topic  has  also  received  considerable  attention,  for  example  the  n-rigid-body 
problem  studied  in  [5],  the  three-rigid-body  problem  discussed  in  [6],  [3],  and  the  effects  of 
finite  extent  of  bodies  on  the  equilibria  of  the  restricted  three-rigid-body  problem  described 
in  [7].  For  the  two- rigid-body  problem  and  the  restricted  two- rigid-body  problem,  there 
are  extensive  discussions  in  [2],  and  in  [23],  [26],  [12],  and  [22]  where  problems  including 
gyrostats  are  considered.  Moreover,  the  steady-state  (rotational)  motions  (or  relative 
equilibria)  of  a  rigid  body  in  a  central  gravity  field  and  their  stability  were  considered 
in  [24],  [21],  [11],  [14],  and  [1].  Similar  problems  are  studied  in  [31]  for  interconnected 
gyrostats.  All  of  the  above  analyses  make  some  approximation  of  the  gravitational  force 
acting  on  the  rigid  body.  We  shall  show  that  for  some  rigid  bodies,  steady-state  motions  of 
the  exact  problem  are  qualitatively,  and  quantitatively  different  from  the  classical  so-called 
Lagrange  or  regular  motions  in  which  one  of  the  principal  axes  coincides  with  the  radius 
vector.  The  precise  differences  are  further  described  below. 

Approximation  of  a  finite  rigid  body  as  a  point  mass  can  be  justified  in  two  distinct 
ways.  A  point  mass  is  an  exact  model  in  the  case  of  a  rigid  body  having  spherical  mass 
distribution  (a  result  due  to  Newton).  Alternatively,  approximation  as  a  point  mass 
appears  to  be  reasonable  if  a  body  is  very  small  in  comparison  to  a  typical  orbit  radius. 
In  the  case  of  an  Earth-satellite  system,  the  justification  for  treating  the  Earth  as  a  point 
mass  is  approximate  spherical  symmetry,  and  the  justification  for  treating  the  satellite  as 
a  point  mass  is  that  it  is  very  small.  Deviations  of  the  Earth  from  spherical  symmetry 
are  an  important  effect  in  low  Earth  orbits,  and  the  problem  has  been  studied  by  many 
workers  including  recent  work  of  Deprit  and  collaborators  [4].  Here  we  are  concerned  with 
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effects  associated  with  the  approximation  of  the  satellite  as  a  point  mass.  Accordingly  we 
shall  hereafter  assume  the  Earth  to  be  exactly  spherical,  and  therefore  exactly  equivalent 
to  (an  extremely  heavy)  point  mass. 

The  equations  of  motion  for  the  restricted  two  body  problem  comprise  a  force  balance 
equation,  and  a  moment  balance  equation.  A  typical  line  of  analysis,  apparently  tracing 
back  to  Lagrange,  and  followed  by  Beletskii  [2],  Sarychev  [25]  and  Roberson  [22]  among 
others,  is  to  retain  the  lowest-order  nonvanishing  term  of  the  potential  in  each  equation. 
As  a  consequence  of  this  strategy  the  force  balance  equation  involves  only  the  zeroth  order 
term  of  the  potential.  Moreover  the  force  balance  equations  decouple  from  the  moment 
balance  equations  to  give  force  conditions  that  are  identical  with  those  of  a  point  mass. 
The  steady-state  solutions,  or  relative  equilibria,  are  circular  orbits  that  are  also  great 
circles  in  the  sense  that  the  center  of  the  circular  orbit  coincides  with  the  center  of  the 
Earth.  Once  the  orbits  have  been  determined,  the  moment  balance  equation  must  still 
be  satisfied.  In  the  moment  balance  the  first  nonvanishing  term  of  the  potential  is  of 
second-order.  It  may  then  be  concluded  that  the  satellite  must  be  oriented  such  that  its 
principal  axes  of  inertia  are  aligned  with  the  radius,  tangent  and  normal  vectors  to  the 
circular  orbit.  When  signs  are  associated  with  the  principal  axes  there  are  (6  X  4  =  )  24 
such  orientations  depending  upon  whether  it  is  the  axis  with  largest  inertia  that  is  tangent 
etc.  Some  of  these  solutions  are  known  to  be  stable,  some  are  known  to  be  unstable,  and 
the  nonlinear  stability  properties  of  others  are  indeterminate. 

In  a  recent  work  [29],  the  above  problem  has  been  reanalyzed  from  the  viewpoint  of 
Hamiltonian  mechanics.  In  particular  it  was  shown  that  the  evolution  equations  can  be 
written  as  a  (noncanonical)  Hamiltonian  system  with  the  potential  being  one  term  in  the 
Hamiltonian.  Different  approximations  of  the  potential  then  provide  different  dynamical 
systems  that  can  be  analyzed.  In  the  case  of  a  second-order  approximation  of  the  potential 
the  classic  results  described  above  were  recovered,  albeit  by  a  different,  more  consistent, 
approximation  scheme.  Interestingly  in  the  analysis  of  the  second-order  model  considerable 
effort  had  to  be  spent  to  exclude  the  possibility  of  non- great- circle  relative  equilibria  in 
which  the  orbit  vector  sweeps  out  a  cone  rather  than  a  disk.  The  possibility  of  non¬ 
great-circle  solutions  is  excluded  a  priori  in  the  approximations  adopted  by  most  previous 
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authors.  In  light  of  this  experience  we  proceeded  to  consider  the  exact  problem  in  which  the 
potential  is  not  approximated,  and  in  this  paper  we  shall  show  that  great-circle  equilibria 
are  exceptional  (i.e.  non-generic)  for  bodies  without  planes  of  symmetry.  Moreover  any 
body  with  no  great-circle  solutions  has  at  least  four  non-great-circle  steady  orbits.  On  the 
other  hand  in  the  exact  problem  for  a  body  with  three  planes  of  symmetry  (as  is  implicitly 
assumed  in  the  second-order  model)  there  are  at  least  24  great-circle  solutions. 

The  counterintuitive  result  concerning  existence  of  non-great-circle  solutions  can 
best  be  explained  by  the  observation  that  as  fax  as  an  inverse  square  gravitational,  field 
is  concerned,  the  center  of  mass  of  a  body  is  not  a  distinguished  point.  However  it  is 
physically  clear  that  the  lateral  displacement  of  the  plane  of  the  orbit  of  the  center  of  mass 
is  necessarily  very  small,  and  in  particular  smaller  than  the  dimension  of  the  satellite. 
Nevertheless  we  will  present  numerical  examples  which  show  that  this  tiny  offset  in  orbit 
can  be  associated  with  a  large  rotation  away  from  the  classic  24  orientations  predicted  by 
the  second-order  model.  Consequently  these  non-great-circle  relative  equilibria  may  be  of 
importance  in  the  understanding  and  design  of  Eaxth  pointing  satellites.  In  particular,  we 
argue  that  large  deviations  in  orientation  are  to  be  expected  when  two  or  more  principal 
inertias  are  nearly  equal. 

In  [2],  [26],  [21]  and  [1],  non-great-circle  relative  equilibria  (sometimes  called  oblique 
regular  motions )  have  already  been  discussed.  Their  approaches  were  based  on  various 
approximate  models  of  the  potential  which  used  local  coordinates  such  as  Euler  angles. 
Here,  we  axe  dealing  with  an  exact  Hamiltonian  model.  The  work  most  closely  related 
is  [1]  in  which  it  was  assumed  that  the  ratio  between  the  body  size  and  the  radius  is 
very  small  and  only  oblique  regular  motions  close  to  the  regular  motions  (i.e.  to  the 
great-circle  relative  equilibria)  axe  observed.  Nevertheless,  we  have  verified  numerically 
that  large  deviations  in  orientation  from  these  classic  motions  are  possible.  In  fact,  we 
exhibit  large  changes  of  orientation  for  a  particularly  simple  idealized  rigid  body  which 
we  call  a  molecule,  namely  six  point  masses,  two  lying  on  each  principal  axis.  In  order 
to  magnify  the  effects  of  higher  order  terms  in  approximation  of  the  potential  function 
(cf.  [14]),  the  three  moments  of  inertia  are  chosen  to  be  close  to  each  other.  On  the 
other  hand,  the  mass  distribution  of  the  molecule  is  designed  to  be  as  asymmetric  as 
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possible  to  emphasize  the  orientation  change.  During  the  numerical  search  for  the  relative 
equilibria  of  such  structures,  continuation  methods  are  adopted  to  solve  efficiently  the 
appropriate  system  of  nonlinear  algebraic  equations  by  successive  application  of  Newton- 
Raphson  iterations.  We  also  apply  an  error  analysis  (due  to  Kantorovich)  which  guarantees 
that  exact  solutions  exist  close  to  our  numerical  approximations.  The  error  analysis  is 
particularly  important  because  the  problems  turn  out  to  be  severely  ill-conditioned  in  the 
appropriate  parameter  regimes.  Indeed  in  order  to  obtain  sufficient  accuracy  to  guarantee 
a  solution,  we  implemented  our  algorithm  in  double-precision  on  a  96-bit  machine  (actually 
in  Fortran  on  a  CRAY  which  leads  to  approximately  29  significant  digits).  The  numerical 
non-great-circle  relative  equilibria  we  found  justify  the  statement  that  while  the  orbital 
offset  is  very  small,  the  change  in  orientation  from  the  regular  motions  can  be  significant. 


2.  Dynamic  Equations  and  Non-dimensionalization 

We  consider  the  motion  of  a  rigid  body  B*  moving  in  a  central  inverse-square 
gravitational  force  field  that  is  generated  by  a  massive  body  Bo* .  The  body  B0*  is  assumed 
to  be  stationary  so  that  a  restricted  two-body  problem  is  at  hand.  A  configuration  of  the 
system  is  depicted  in  Figure  2.1.  We  assume  that  the  reference  or  inertial  frame  is  located 
at  the  center  O  of  the  inverse  square  field  due  to  the  mass  of  Bo* .  Let  C  denote  the 
center  of  mass  of  the  moving  body  B* ,  and  let  r*  denote  the  vector  from  0  to  C  in 
the  inertial  frame.  The  attitude  of  the  rigid  body  B*  is  described  by  a  set  of  orthogonal 
basis  vectors  located  at  C  and  fixed  relative  to  B* ,  which  will  be  termed  the  body  frame  . 
Let  the  orientation  transformation  from  the  body  frame  to  the  inertial  frame  be  denoted 
by  B ,  which  is  a  rotation  matrix  or  element  of  the  special  orthogonal  group  50(3).  By 
elementary  analysis,  the  time  derivative  of  the  rotation  matrix  can  be  written  as 

B  =  BQ*,  (2.1) 

where  ~  maps  the  instantaneous  angular  velocity  0*  €  IR3  of  B  (or  B*)  in  the  body 
frame,  to  an  element  in  so(3),  the  space  of  3  x  3  skew-symmetric  matrices,  according  to 
the  rule 
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(2.2) 


/or i\  /  o  -ft*3  ft*2  \ 

fi*2  =  ft*  3  0  -ft*i 

\ft*3/  \-ft*2  0*1  0  / 


Figure  2.1.  A  Rigid  Body  in  a  Central  Gravitational  Field 


In  (2.1)  and  what  follows,  we  adopt  the  convention  that  variables  with  an  asterisk  * 
have  dimension,  and  the  dot  operator  '  denotes  the  derivative  of  the  corresponding  variable 
with  respect  to  the  dimensional  time  variable  t* . 

Let  Q*  be  the  vector  from  C  to  an  arbitrary  material  point  in  B*  relative  to  the  body 
frame,  and  denote  the  associated  mass  element  by  dm(Q*).  Then  the  linear  momentum 
of  B*  can  be  written  as 

P*  =  [  4 (r*  +  BQ *)  dm(Q*)  =  (2.3) 

J  B*  dt 

where  m  is  the  total  mass  of  B* .  Here  the  fact  that  fg,  Q*dm(Q*)  =  0  has  been  used.  On 
the  other  hand,  the  angular  momentum  of  B*  about  0  (expressed  in  the  inertial  frame)  is 

=  f  (r*  +  BQ')  x  -Lr*  +  BQ')  dm(Q') 

Jb*  dt  (2.4) 

=  r*  x  mr*  +  HI*  ft*, 

where 
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Q*Q*  dm(Q*) 


(2.5) 


is  the  moment  of  inertia  of  B*  relative  to  the  body  frame. 

The  gravitational  force  exerted  on  B*  by  Bo*  can  be  derived  from  the  potential 

G 


Ibq *  Ib 


b*  k*  ~  (r*  +  BQ*)\ 


dm(Q*)dmo(q*), 


(2.6) 


where  G  is  the  universal  gravitational  constant,  and  q*  is  a  vector  from  O  to  a  material 
point  in  Bo*  with  mass  dm0(q*).  Now  we  assume  that  the  body  Bo*  has  spherical 
symmetry.  Then  by  elementary  analysis  it  can  be  shown  that  the  potential  function  V 
simplifies  to 

GM 


V 


=  -L 


&  jr*  +  BQ* 


dm(Q*), 


(2.7) 


where  M  is  the  total  mass  of  the  body  Bo* .  Consequently,  assuming  there  are  no  torces 
external  to  the  two-body  system,  the  resultant  force  on  B*  can  be  expressed  as 

GM(r*  +  BQ*) 


f 


~J. 


*  Ir*  +  BQ*  |* 


dm(Q*). 


(2.8) 


Similarly,  assuming  there  is  no  external  couple  acting  on  B* ,  the  resultant  torque  is 


f  /  *  GM l 

T  =  -l(r  +BQ-)x  — 


GM(r*  +  BQ*) 


+  BQ 


*13 


dm(Q*)  =  0. 


(2.9) 


Balance  of  linear  momentum  then  implies,  from  (2.3),  and  (2.8), 


f  GM\ 


GM{r*  +  BQ*) 


+  BQ*  |3 


dm(Q*). 


(2.10) 


On  the  other  hand,  balance  of  angular  momentum  implies 


7T*  =  0. 


(2.11) 


Let  II*  =  1*0*  denote  the  instantaneous  body  angular  momentum  of  B* .  Then  from  (2.1), 
(2.4),  and  (2.9),  we  can  derive  an  evolution  equation  for  II* 

.-IT,.  ,  r  GM(BTr-  X  CT) 


n-  =  n'  xr‘ir  + 


L 


s-  k*  +  BQ*  |3 


dm(Q*). 


(2.12) 
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Consequently,  equations  (2.1),  (2.3),  (2.10),  (2.12)  completely  describe  the  motion  of  the 
rigid  body  B*  moving  in  the  gravitational  force  field. 

Now  we  define  convected  or  body  variables, 

A*  =  BTr\  n*  =  BTp\ 


which  axe  vectors  expressed  in  coordinates  with  respect  to  the  moving  body  frame.  In 
terms  of  these  variables,  the  equations  of  motion  may  be  written  as 


‘xn*  +  f 

JB 


n*  =  n*  x  r 

A*  =  A*  x  I*-1  IT  +  p*/m, 

/> 


GM(\*  X  or) 
b*  |A*  +  Q*  I3 


B  = 


/**  x  r-'ir  -  / 

JB 


GM(\*  +  Q*) 
5-  I  a*  +  Q* I3 


dm(Q*), 


dm(Q*), 


(2.13a) 

(2.136) 

(2.13c) 

(2.13d) 


We  now  nondimensionalize  this  system.  With  mass,  length,  and  time  scales, 


m. 


l  = 


'fr(I*) 


m 


and  T 


l3 


V 

GMJ  ’ 


the  nondimensional  variables  axe 


,  A* 

A  =  — ,  n 


mlT~1  ’ 


„  n*  '  ,  r 

n  =  “d  *  =  f • • 


(2.14) 


(2.15) 


Next  we  introduce  a  nondimensional  body  (or  domain  of  integration)  #.  From  the 
dimensional  mass  measure  dm ,  we  define 

(2.16) 


Q  =  MQ)  = 

j  .  m 


and  accordingly, 


[  dv{Q)  =  /  -*n(Q*)  = 
XB  ./b*  m 


(2.17) 


The  nondimensional  moment  of  inertia  of  B  is 


I  = 


fr(I*) ’ 


(2.18) 
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so  that 


tr(  I)  =  1.  .  (2.19) 

In  terms  of  these  nondimensional  variables,  the  dynamical  equations  (2.13)  can  be 
expressed  as 


n'  =  n  x  r’n  +  J  A  *  ®  *-(<?),  (2.200) 

A'  =  Ax  I^H  +  /i,  (2.206) 

/*'  =  p  x  i_,n  -  J  |2-±A  &4Q),  (2.20c) 

S'  =  (2.20d) 


when  the  prime  '  denotes  a  derivative  with  respect  to  the  nondimensional  time  t .  It  is 
the  nondimensional  system  (2.20)  that  will  be  investigated  in  the  following  Sections. 


3.  Relative  Equilibria 

Equations  (2.20  a,b,c)  do  not  include  the  attitude  B  and  are  therefore  decoupled  from 
(2.20d).  Accordingly,  the  equations  of  motion  (2.20  a,b,c)  on  the  nine  dimensional  space 
(II,  A,  fi)  €  IR9  will  be  called  the  reduced  dynamics.  Equilibria  of  the  reduced  dynamics 
give  rise  to  motions  in  the  original  space  (r*,  B,p*,  11*)  satisfying  (2.1),  (2.3),  (2.10),  (2.12) 
such  that  the  rigid  body  B *  exhibits  a  steady  spin  about  the  fixed  center  O,  while  the 
center  of  mass  C  moves  in  a  circular  orbit.  Such  a  motion  will  be  referred  to  as  a  relative 
equilibrium.  All  such  motions  are  equilibria  relative  to  a  steadily  rotating  frame.  There  is 
no  a  priori  reason  that  the  center  of  the  circular  orbit  traced  by  C  need  coincide  with  the 
origin  O .  Indeed  we  will  exhibit  relative  equilibria  where  O  is  not  the  center  of  the  orbit. 

Conditions  of  relative  equilibria  are  easily  found  from  the  reduced  dynamics.  In 
terms  of  the  reduced  variables,  we  have  the  following  equations, 

nxl~'n  +  JsWTwMQ)  =  0’  (31o) 
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(3.14) 

(3.1c) 


a  x  i  1n  +  fi  =  o, 

<ixI'ln  -  Lwrw**®  =  °- 


Discussion  of  the  existence  of  solutions  to  (3.1),  and  their  qualitative  and  quantitative 
features  axe  the  main  theme  of  this  paper. 

In  [29],  a  Hamiltonian  approach  to  this  problem  was  exploited.  It  was  shown  that 
the  reduced  dynamics  (2.20  a,  b,  c)  can  be  written  as  a  Hamiltonian  system 

4  |  A  I  =  AVI.  (3.2) 


with  Hamiltonian 

jfCn.A.rt  =  i  <  n,i-n  >  ^  *tQ).  (**) 

and  Poisson  tensor, 

/n  a  A 

A  =  A  0  1  ,  (3.4) 

\£  -1  0  J 

which  is  a  9  x  9  skew-symmetric  matrix  with  block  structure.  Here  the  notation  is 
defined  through  (2.2),  and  1  denotes  the  3x3  identity  matrix.  The  Poisson  tensor  A  is 
rank-degenerate.  In  fact,  its  null  space  is  one- dimensional  and  is  spanned  by  the  vector 

/  n  +  Ax/i  \ 

VC(n,A,/x)  =  Ux(H  +  Ax/i)  ,  (3.5) 

\(U  +  X  x  n)  x  X  J 


where  C  is  the  total  angular  momentum, 

c(iu,p)  =  lin  +  Ax^p. 


(3.6) 


The  function  C  is  sometimes  called  a  Casimir  function  of  the  noncanonical  Hamiltonian 
structure,  cf.  [9].  It  is  clear  from  (3.2)  that  at  relative  equilibria,  the  vector  VH  must  lie 
in  the  null  space  of  A,  in  other  words, 

Vff( n,A,/i)  =  C  VC(n,A,^),  (3.7) 
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for  some  scalar  (.  Accordingly,  there  is  a  variational  principle  which  characterizes  the 
relative  equilibria,  namely 

make  stationary  H( II,  A,  fi) 

(3.8) 

subject  to  the  constraint  C(II,  A,  fx)  =  constant. 

The  first-order  conditions  for  (3.8)  are  nothing  other  than  (3.7),  which  can  also  be  expressed 


as 


/  i-!n  \  /  n  +  A  xix 

vAy(A)  =  c  U  x  (n  +  a  x  n) 

\  H  )  \(Tl  +  X  x  fx)  x  \ 


(3.9) 


where 


y(A)  =  -  f 
Jb 


b  |A  +  Q\ 


dm(Q). 


(3.10) 


Equation  (3.9)  can  be  immediately  simplified  to 


('?’) 


C(n  +  a  x  fx) 
(x  xi-'n 
i-Jn  x  a 


(3.11) 


Thus  ix  can  be  eliminated  from  (3.11)  to  yield 

(  i-xn  \  _  f  c(n  +  a  x  (i_1n  x  A)) 
VvAy(A);  “  ^  (i-1n  x  A)  x i_1n 

By  introducing  the  body  angular  velocity  ft  =  I-1  II,  we  find 

Ift  +  A  x  (ft  x  A)  =  -ft, 


(3.12) 


(3.13) 


(1)x1)x8-VjK(A)  =  0. 

It  is  easy  to  check  that  (3.13)  are  the  first-order  conditions  for  the  variational  principle, 

(3.14) 


where 


make  stationary  17(ft,  A)  subject  to  —  |ft|2  =  c, 

Jit 


U  =  \  <  ft,  Ift  >  +  \  |ft  X  A|2  -  y(A). 


(3.15) 


It  is  this  variational  principle  that  provides  the  starting  point  of  our  proof  of  the  existence 
of  relative  equilibria.  By  further  rearrangement,  (3.13)  can  be  rewritten  as 
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(i  +  xTx)o,  =  pn, 

(3.16a) 

Ax  ~  L  ia  +  <?p ivm 

(3.166) 

5  lap  =  c. 

(3.16c) 

where  ATA  and  ClTCl  are  3x3  matrices  constructed  from  the  vectors  A  and  fl  through 
formula  (2.2).  By  construction,  solutions  of  (3.1)  and  solutions  of  (3.16)  are  in  one-to-one 
correspondence. 

It  should  be  remarked  that  the  Hamiltonian  system  (3.2)  is  the  Poisson  reduced 
system  of  (2.20)  in  the  sense  of  Poisson  reduction  [13].  Relative  equilibria  characterized 
by  (3.1)  or  (3.16)  exactly  coincide  with  the  notion  of  relative  equilibrium  associated  with 
the  Poisson  reduction.  Accordingly,  the  principle  of  symmetric  criticality  [19]  [2,8]  can  be 
applied,  namely,  the  configuration  components  of  the  relative  equilibria  are  the  critical 
points  of  the  augmented  potential  function,  which  for  our  problem  is  —  U  .  Consequently, 
the  variational  principle  (3.14)  is  a  particular  manifestation  of  the  principle  of  symmetric 
criticality. 


4.  Great-circle  Relative  Equilibria 

If  we  make  the  further  assumption  that  the  rigid  body  B  has  spherical  symmetry, 
then  it  can  be  checked  that  the  dynamical  behavior  of  the  center  of  mass  C  is  exactly 
that  of  the  Keplerian  point  mass  model.  It  is  well  known  that  for  a  point  mass  moving  in 
a  central  force  field,  the  motion  is  constrained  to  a  plane  passing  through  the  origin  O . 
However,  for  an  arbitrary  rigid  body,  this  may  not  be  the  case.  In  particular,  at  relative 
equilibrium,  the  radius  vector  from  the  origin  O  to  the  center  of  mass  C  may  either  trace 
a  disc,  or  trace  a  cone,  cf.  Figure  4.1(a),  (b).  Those  relative  equilibria  in  which  C  traces 
a  circle  that  resides  on  a  plane  passing  the  origin  O  will  be  termed  great- circle  relative 
equilibria.  Otherwise,  they  will  be  called  non- great- circle  relative  equilibria.  In  this  section, 
we  study  only  the  great-circle  relative  equilibria. 

From  Figure  4.1(a),  it  is  clear  that  the  necessary  and  sufficient  condition  for  great- 
circle  relative  equilibria  is 
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vj3 

I 

I 

I 


(b) 

Figure  4-1-  (a)  Great-circle  vs.  (b)  Non-great-circle  Relative  Equilibria. 

Here  r  =  BX  is  the  radius  vector  in  inertial  frame,  and  £  =  B£l  'points 
along  the  axis  of  rotation.  The  classic  approximate  analysis  predicts  great- 
circle  orbits  in  which  the  orbit  vector  r  sweeps  out  a  dish.  However,  the 
Hamiltonian  approach  to  the  exact  problem  for  asymmetric  bodies  predicts 
non-great-circle  solutions  in  which  r  sweeps  out  a  cone. 

A- H  =  0,  (4.1) 

namely,  the  axis  of  rotation  is  perpendicular  to  the  radius  vector  of  C .  Notice  that  equation 
(3.18a)  can  be  rewritten  as 

(I  +  | A j2 1  -  AAr)ft  =  0  SI,  (4.2) 
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where  1  denotes  the  (3x3)  identity  matrix  and  \\T  is  a  (3  x  3)  rank  one  matrix.  Thus 
for  a  great-circle  relative  equilibrium,  in  which  A  •  Q  =  0 ,  we  must  have 

in  =  (/?  -  |a|2)  n,  (4.3) 


or  equivalently  Cl  must  be  an  eigenvector  of  the  moment  of  inertia  tensor.  Moreover, 
condition  (3.16b)  reduces  to 


\Cl\2X 


r  o  +  0) 

B  |A  +  Ql3 


dm(Q). 


(4.4) 


By  regarding  |fi|2  in  (4.4)  as  a  Lagrange  multiplier,  it  can  be  shown  that  (4.4)  is  the 
first-order  condition  for  the  variational  problem, 


make  stationary  V(A)  subject  to 


2  'A|2  -  S’ 


(4.5) 


where  V  is  defined  in  (3.10).  In  summary,  in  order  for  (Q,  A,/z)  to  be  a  great-circle 
relative  equilibrium,  the  vector  A  must  be  a  solution  of  the  problem  (4.5),  and  Cl  must  be 
an  eigenvector  of  the  inertia  tensor,  i.e.  a  principal  axis.  With  condition  (4.1),  we  then  see 
that  A  must  also  lie  in  a  principal  plane  (i.e.  a  plane  perpendicular  to  a  principal  axis). 

The  argument  can  be  reversed.  If  there  is  a  vector  Ao  solving  the  variational  problem 
(4.5),  there  exists  a  multiplier  k  such  that 

VaV-(Ao)  =  « Ao.  (4.6) 


Moreover,  for  physically  relevant  situations  the  rigid  body  does  not  overlap  the  center  of 
the  field.  Thus 

|A|2  +  A  -Q  >  |A|2  -  |A||Q|  =  |A|(|A|-|Q|)  >  0.  (4.7) 

It  can  then  be  proved  from  (4.6)  and  (4.7)  that  the  multiplier  k  is  necessarily  positive. 
By  comparing  (4.6)  and  (3.16b)  and  if  A  lies  in  a  principal  plane,  we  can  construct  a 
great-circle  solution  by  picking  |Q|2  =  «,  and  Cl  the  principal  axis  with  Cl  ■  A  =  0.  In  fact, 
there  are  at  least  two  such  solutions  corresponding  to  the  choices  Cl  and  —Cl .  Accordingly, 
if  (4.6)  holds  for  a  A0  in  a  principal  plane,  we  can  find  at  least  two  great-circle  relative 
equilibria.  Furthermore,  any  plane  of  symmetry  of  the  body  B  is  necessarily  a  principal 
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plane,  and  there  are  at  least  two  critical  points  of  F(A)  in  the  symmetry  plane  (a  maximum 
and  a  minimum).  As  a  consequence,  we  have  the  following  theorem. 

THEOREM  4.1. 


For  a  rigid  body  having  a  plane  of  symmetry,  there  are  at  least  four  great- circle 
relative  equilibria.  Furthermore,  if  the  rigid  body  is  symmetric  with  respect  to  two  planes, 
there  are  at  least  eight  great-circle  relative  equilibria,  and  for  a  rigid  body  with  three 
planes  of  symmetry,  there  are  at  least  24  great-circle  relative  equilibria. 

I 

The  second  and  third  statements  in  Theorem  4.1  follow  from  similar  arguments  as 
those  described  above.  The  classic  second-order  approximate  models  analyzed  in  [n]  im¬ 
plicitly  assume  that  the  body  has  three  planes  of  symmetry,  and  they  find  exactly  24 
grea-  -circle  relative  equilibria. 


5.  Non-great-circle  Relative  Equilibria 

In  this  section,  we  shall  establish  the  existence  of  non-great-circle  relative  equilibria 
for  non-symmetric  bodies.  We  observed  in  the  previous  section  that  at  a  great-circle 
relative  equilibrium,  the  variational  problem  (4.5)  must  have  a  solution  A  lying  in  a 
principal  plane.  This  fact  can  be  used  to  prove  that  for  some  particular  V  (i.e.  for 
some  particular  rigid  body  B),  there  are  no  great-circle  relative  equilibria,  cf.  [29],  [27]. 
In  particular,  let  the  rigid  body  be  a  “molecule”  consisting  of  six  point  masses,  two  on 
each  principal  axis,  cf.  Figure  5.1.  (The  molecule  with  unequal  principal  moments  of 
inertia  is  an  example  in  Beletskii  [2]  satisfying  the  conditions  for  the  existence  of  relative 
equilibrium.)  It  was  shown  in  [27]  that  if  x\  ^  xi ,  y\  ^  t/2  >  and  z\  ^  z<i ,  there  is  no  great- 
circle  relative  equilibrium.  Now  if  we  could  prove  the  existence  of  solutions  to  the  problem 
(3.14),  then  for  the  cases  in  which  there  is  no  great-circle  relative  equilibrium,  there  must 
be  a  non- great- circle  relative  equilibrium.  In  fact,  we  shall  show  that  for  a  general  body, 
variational  principle  (3.14)  is  not  useful  for  proving  existence  of  relative  equilibria.  The 
difficulty  is  associated  with  “critical  points  at  infinity” ,  namely  solutions  in  which  $7  and 
A  are  parallel  with  |A|  — +  00.  Fortunately,  this  difficulty  can  be  finessed  due  to  the  pure 
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quadratic  dependence  upon  fi.  We  introduce  another  variational  principle,  namely 

make  stationary  U(Q.,  A)  =  -  <  Ifl  >  +  -  \Q.  x  A|2, 

2  2 


subject  to 


i  |0|2  =  ci,  and  y(A)  =  c2, 


(5.1) 


where  we  have  introduced  an  additional  artificial  constraint  with  c2  lying  in  the  interval 
(--  jB(l/\Q\)dv{Q),  0).  The  Lagrangian  associated  with  (5.1)  is 


l0  =  i  <  n,  in  >  +  i  |n  x  A|>  -  |np  -  c.)  -  «v(A)  -  c2).  (5.2) 


m3 


Figure  5.1.  Molecule 


Solutions  of  (5.1)  with  /?2  =  1  correspond  to  solutions  of  (3.14)  for  some  value  of  c  and 
therefore  to  relative  equilibria.  The  first-order  conditions  associated  with  (5.2)  are 


(i  +  xTX)n  =  Pi  a, 

(5.3  a) 

CiTCix  -  p2VxV{\)  =  o, 

(5.36) 

=<=„ 

(5.3c) 

P 

II 

3 

(5.3d) 

By  taking  the  scalar  product  of  A  with  (5.3b),  we  get 
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(5.4) 


ft  /B(|A|A  '+'Q|^)  dm(Q)  =  (MW  -  (n'A>2)- 

Notice  that  the  right  hand  side  of  (5.4)  is  always  nonnegative.  On  the  other  hand,  from 
(4.7),  it  is  seen  that  the  multiplier  02  of  solutions  of  problem  (5.1)  must  be  nonnegative, 

02  >  0,  (5.5) 


with  the  equality  holding  if  and  only  if  ft  and  A  axe  parallel. 

We  now  compare  conditions  (5.3)  and  (3.16).  It  is  clear  that  if  there  is  a  solution 
(ft,  A,  0i,  02)  of  (5.3)  with  02  =  1,  then  (3.16)  is  solved  by  (ft,  A,  0\)  with  c  ~  c\. 
Moreover,  since  (5.3a)  is  linear  in  ft,  we  can  scale  ft  arbitrarily  without  violating  (5.3a). 
Accordingly,  if  we  have  a  solution  of  (5.3)  such  that  02  >  0,  we  obtain  a  solution 
(ft/Vfo)  A,  0i)  of  (3.16)  with  c  =  C\j02.  On  the  other  hand,  the  zeroth-order 
approximation  of  (5.3b)  implies  that  02  ~  2c\/c\  (Kepler’s  frequency  formula).  Thus 
by  holding  Ci  fixed,  and  by  varying  C2 ,  we  can  achieve  arbitrary  02 .  The  constant  c 
is  then  also  seen  to  be  arbitrary.  As  a  consequence,  in  order  to  prove  the  existence  of 
solutions  for  problem  (3.14)  for  some  c,  it  suffices  to  find  a  solution  of  the  problem  (5.1) 
with  02  >  0 . 

Next,  we  consider  those  extremals  at  which  02  —  0.  From  (5.5),  we  see  that  for  such 
extremals  ft  is  parallel  to  A ,  which,  in  turn,  shows  that  ft  must  be  an  eigenvector  of  I , 
say 

Ift  =  lift, 


or  0i  =  Ii ,  cf.  (5.3a).  It  is  then  easy  to  check  that  the  objective  function  [7(ft,  A)  cannot 
assume  its  minimum  value  at  such  extremal  points.  In  fact,  the  Hessian  of  Lg  for  02  =  Q 
can  be  found  to  be 


D(Cl,\)LU 


I  +  ArA-Iil  2Aft  —  ftA  \ 
2ftA  -  Aft  ftrft  )  ’ 


(5.6) 


and  the  tangent  space  at  (ft.  A)  to  the  constrained  subspace  can  be  expressed  as, 


v  =  {(yi,y2)elR6  :  ft  •  yi  =  0,  VaK(A)  •  y2  =  0  }. 
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Since  1)  is  a  principal  axis  of  the  rigid  body  B ,  the  vector  y\  must  be  in  a  principal  plane. 
We  compute 

(yf  °)  Dt„ML0(n,X) 

=  yiTIyi  +  (A  x  yi)2  -  Ji|yi|2 
>  (|A|2  -  (Ii  -min{J,-}))  jyx|2, 

since  we  also  have  A  •  yi  =0.  Moreover,  through  the  non-dimensionalization  process,  we 
have  tr(I)  =  1,  which  implies  \Ii  —  min{J,}|  <  1.  Thus  with  the  natural  assumption 
|A|  >  1,  cf.  (4.7),  the  difference  |A|2  —  (7i  —  min{/j})  is  always  positive.  Consequently, 
the  quadratic  form  of  the  Hessian  restricted  to  V  can  never  be  negative-semidefinite.  This 
implies  that  the  critical  point  with  /?2  =  0  can  never  be  a  maximizer  for  the  problem 
(5.1).  However,  problem  (5.1)  has  a  continuous  objective  function  U  on  a  closed  and 
bounded  constraint  set.  It  follows  that  there  must  be  a  maximizer  for  the  problem..  The 
above  argument  reveals  that  this  maximizer  must  be  a  critical  point  with  fa  >  0.  As  a 
consequence,  there  must  be  a  solution  for  the  problem  (5.4).  Recall  that  for  the  molecule 
example,  cf.  Figure  5.1,  there  is  no  great-circle  relative  equilibrium.  Thus  the  solution 
whose  existence  is  proven  here  must  give  rise  to  a  non-great-circle  relative  equilibrium. 

More  detailed  information  can  be  obtained  by  applying  Morse  theory  to  this  problem. 
The  variational  principle  (5.1)  can  be  re-stated  as  finding  critical  points  of  U  on  the 
manifold  HP2  x  S2 ,  by  modding  out  the  ZS2  symmetry  on  fi.  It  can  be  checked  that  the 
Betti  numbers  for  the  manifold  HP2  x  S 2  are  /?0  =  1 ,  /?i  =  1 ,  02  =  2 ,  0z  =  1 ,  04  =  1 . 
On  the  other  hand,  the  indices  of  the  critical  points  with  02  =  0  can  be  found  by  counting 
the  number  of  negative  eigenvalues  of  the  Hessian  (5.6)  restricted  to  the  tangent  space  V . 
It  turns  out  that  there  are  two  solutions  for  each  index  0,  1,  and  2.  Now  invoking  the 
Morse  inequalities  [16],  we  conclude  there  is  at  least  one  solution  for  each  of  the  indices 
3  and  4,  respectively.  Consequently,  by  including  the  2£2  symmetry,  it  may  be  concluded 
that  there  are  at  least  four  non-great-circle  relative  equilibria. 

In  [1],  a  proof  for  the  existence  of  non- great-circle  relative  equilibria  for  approximate 
models  was  given  based  on  the  assumptions  that  the  principal  central  moments  of  inertia 
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are  not  equal,  and  the  body  size  is  sufficiently  small  in  comparison  with  the  distance  to 
the  center  of  the  field.  The  variational  proof  given  above  applies  to  the  exact  model  as 
well  as  to  its  approximations  and  makes  less  restrictive  hypotheses. 


6.  A  Numerical  Search  and  Error  Analysis 

The  existence  of  non-great- circle  relative  equilibria  for  certain  bodies  B  are  ensured 
by  the  arguments  presented  in  the  previous  section.  We  now  seek  such  relative  equilibria 
numerically,  and  justify  the  statement  that  while  the  offset  of  the  circular  orbit  is 
necessarily  very  small,  the  associated  attitude  change  may  be  quite  significant.  The 
numerical  method  adopted  and  the  associated  error  analysis  are  described  in  this  Section. 
Specific  examples  will  be  given  in  Section  7. 

We  consider  the  example  of  a  molecule  as  depicted  in  Figure  5.1.  Recall  that  the 
conditions  for  relative  equilibria  are  given  in  (3.16),  which  can  be  rewritten  as 

g 

ini’*  -  (Q-A)n  -  ^E"[x+^.r'  =  °’  (6-la) 

K!  +  |A|2fi  -  (Sl-A)A  +  =  0,  (6.16) 

1  -  c  =  0,  (6.1c) 

where 

>  Qz 

Qa  =  £4^  ,  Q5  —  ^  0  ^  ,  Qe 

and  A  is  a  parameter  included  for  numerical  purposes.  Physically  the  parameter  A 
is  associated  with  scaling  in  time.  It  will  be  chosen  later  to  minimize  numerical  ill- 
conditioning.  Equation  (6.1)  comprise  a  nonlinear  system  of  seven  equations  in  seven 
unknowns  (A ,  0,/3)  G  JR7.  A  numerical  continuation  method  will  be  used  to  find  solutions 
of  these  equations,  and  therefore  to  yield  relative  equilibria. 
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First  we  make  some  observations.  It  has  been  proven  that  if  a  molecule  is  asymmetric, 
there  is  no  great-circle  relative  equilibrium.  Thus  for  the  purpose  of  illustrating  examples 
of  non-great-circle  relative  equilibria  it  suffices  to  consider  bodies  of  this  very  special  type. 
On  the  other  hand,  if  the  molecule  is  symmetric  (i.e.  x\  =  x2  >  x3  =  XA »  £5  =  ),  there 

are  24  great-circle  relative  equilibria  which  are  known  explicitly.  Such  explicit  solutions 
will  be  the  starting  points  for  our  continuation  methods.  Consider  the  following  case 


fl  = 
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because  of  the  symmetry  of  the  molecule.  Thus  (6.1a)  can  be  reduced  to 


ab 2  —  A(  - 


mi 


+ 


m2 


+ 


2rri3a 


+ 


2m5<z 


(a  +  x  1)2  (a -2:2)2 

Similarly  for  solutions  of  the  form  (6.2),  equation  (6.1b)  reduces  to 

+  a2b  —  j3b  =  0. 


:)  =  0, 


(6.2) 


(6.3) 


With  a,  6  =  \/2c,  and  /?  satisfying  (6.3)  and 

/?  =  -(J2  +  a2),  (6.4) 

then  A,  0  of  the  form  (6.2)  provide  a  solution  of  (6.1),  and  are  thus  a  relative  equilibrium. 
Similar  arguments  can  be  applied  to  ail  other  great-circle  solutions.  With  the  solutions 
for  symmetric  cases  as  starting  points,  it  is  possible  to  follow  the  path  of  solutions  by 
numerical  continuation  as  the  parameters  are  varied  to  the  desired  asymmetric  case. 

The  numerical  continuation  or  path-following  method  is  discussed  in  many  places, 
see  e.g.  [loj  [30]  [8],  and  references  therein.  The  basic  idea  is  to  solve  a  sequence  of 
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problems  as  a  parameter  is  varied,  starting  with  a  known  solution,  so  that  at  each  step  a 
locally  convergent  algorithm  can  be  applied  to  yield  the  new  solution,  which  will  form  the 
starting  point  for  the  next  iteration.  A  homotopy  provides  the  linkage  between  the  known 
simple  solutions,  and  the  desired  answer.  Since  in  our  problem  the  relative  equilibria 
for  a  symmetric  molecule  are  simply  obtained  by  solving  (6.3),  it  is  natural  to  start  the 
homotopy  from  these  symmetric  cases.  Let  (6.1)  be  represented  as 


F(mi,m2,m3,m4,m5,m6,  A,f2,/3)  =  0. 


(6.5) 


For  a  completely  asymmetric  molecule,  we  have  mi  ^  m2,  m3  ^  m4 ,  and  m5  7^  m6 . 
Define 


m,  = 


mi  +  m2 


~2 - ’  = 


m3  +  m4 


m5  +  m6 


Three  homotopies  will  be -used  in  series  to  fulfill  the  path-following  process,  namely, 

JTi(A,ft,/?,r)  =  F(rmi  +  (1  -  r)mI,rm2  +  (1  -  r)  mx,  my,  my,  mz,  mz,  A,  ft,  /?),  (6.6a) 
£T2(A,D,/3,r)  =  F(m1,m2,rm3  +  (1  -  r)mJ/,rm4  +  (1  -  T)my,mz,mg,\,$l,0),  (6.66) 
Ff3(A,0,/3,T)  =  F(m1,m2,m3,m4,rm5  +  (1  —  r)mz,rm6  +  (1  —  r)mz,  A ,fi,/?).  (6.6c) 


By  slowly  varying  r  at  each  stage,  locally  convergent  algorithms  such  as  the  Newton- 
Raphson  method  lead  to  the  desired  solutions.  In  order  to  successfully  find  the  solutions 
of  (6.3),  we  note  the  fact  that  for  A  •  Q.  =  0  and  |A|  >  max{x;,  i  =  1,  •  •  • ,  6} , 


|A|W 


A, 


(6.7) 


which  is  the  Kepler’s  frequency  formula  for  the  zeroth  order  approximate  model,  cf.  [29]. 
Accordingly,  for  given  A,  (6.3)  can  be  solved  by  choosing  the  initial  value  (A/62)1/'5 . 

The  above  path-following  method  will  later  be  termed  method  of  continuation  in 
mass.  We  now  study  the  local  numerical  algorithm  for  finding  solutions  at  each  step  along 
the  path.  Due  to  the  smallness  of  the  ratio  between  the  body  size  (the  X{ )  and  the  radius 
j  A | ,  the  numerical  search  for  solutions  of  (6.1)  becomes  highly  ill-conditioned.  The  constant 
A  plays  an  important  role  in  making  our  numerical  methods  tractable.  Before  performing 
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the  error  analysis,  we  describe  a  theorem  of  Kantorovich  which  will  play  a  pivotal  role  in 
proving  that  actual  solutions  lie  close  to  our  numerical  approximates. 

THEOREM  6.1.  (  Kantorovich) 

Let  D  C  1R"  be  an  open  set,  F  :  D  — ►  1R”  be  continuously  differentiable,  and 
the  matrix  F'  be  Lipschitz  continuous  in  D  with  Lipschitz  constant  a.  Let  xo  G  D , 
with  F'(xq)  nonsingular,  and  suppose  there  are  positive  constants  a  and  7  such  that 
||F'(x0)-1||  <  cr,  ||i?'(i0)-1  jP(aro)H  <  7,  and  h  =  aa-y  <  1/2.  Set  e*  =  (1— vT"—  2 h)/aa. 
Let  B(xo>e*)  be  the  open  ball  with  center  xo  and  radius  e*,  and  suppose  further  that  h 
is  so  small  that  J5(xo,e*)  C  D.  Then  there  is  a  unique  x*  €  B(x o,e*)  with  F(x*)  =  0. 

I 

The  proof  of  this  theorem  can  be  found  in  [10],  [18],  or  [17].  The  first  part  of  the  theorem 
provides  a  way  to  measure  the  error  bound  for  the  approximate  solutions,  and  could  also 
serve  as  a  criterion  for  convergence.  In  particular,  if  an  approximate  solution  satisfies  the 
conditions 

7<e,  CN77  <  1/2,  (6.8) 

then  the  error  of  this  approximate  solution  is  bounded  by  e,  As  a  consequence, 
condition  (6.8)  can  be  used  to  guarantee  the  convergence  of  an  algorithm. 

In  order  to  apply  Theorem  6.1  to  our  problem,  a  value  of  A  (appearing  in  (6.1)) 
muct  be  chosen  and  there  are  some  constants  to  be  determined.  The  detailed  derivation 
of  these  constants  is  presented  in  Appendix  B.  Here  we  summarize  the  results  by  noting 
that  the  conditions  on  A  are  given  by,  cf.  (B.9),  (B.15), 

A  ~  |A|4,.  and  A  <  -J ^ — ,  (6.9) 

*  ^max 

and  suitable  forms  for  the  constants  appearing  in  Theorem  6.1  are,  cf.  (B.3),  (B.12), 
and  (B.17), 

Yl 

7  =  max{&ri,  i  =  l,  a  =  - ,  and  a  ~  2n|A|.  (6.10) 

Pmin 

Here  6x  solves  (B.ll),  and  pmin  is  the  minimum  eigenvalue  of  Fr(\,Sl,  /3) .  With  those 
constants  determined  through  (6.10),  criteria  (6.8)  can  be  applied  to  the  Newton- Ratphson 
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method  to  find  approximate  solutions  with  an  associated  error  bound.  The  constant  A 
determined  from  (6.9)  makes  it  possible  to  deal  with  numerical  ill-conditioning.  These 
considerations  will  be  further  described  in  the  next  section. 


7.  Examples 

In  this  section,  the  ideas  described  in  Section  6  are  combined  to  generate  an  algorithm 
which  finds  non-great-circle  relative  equilibria.  Here,  we  only  look  at  simple  models,  namely 
a  rigid  body  B  with  the  structure  of  a  molecule,  but  in  principle,  the  methods  we  use  could 
be  extended  to  more  complicated  mass  distributions.  In  the  first  example,  the  parameters 
of  the  molecule  are  motivated  by  the  parameters  of  the  motion  of  the  satellite  Phobos 
around  the  planet  Mars.  Non-great-circle  relative  equilibria  with  small  deviations  from 
the  classic  solutions  axe  axe  numerically  found.  We  then  fix  the  size  of  the  molecule  and 
the  radius  of  the  orbit  |A|,  and  study  the  effects  of  different  moments  of  inertia  on  the 
relative  equilibria.  The  second  example  is  artificial,  but  does  include  the  parameter  regimes 
in  which  a  space  station  orbiting  the  Earth  would  fall,  and  serves  as  a  demonstration 
of  large  orientation  change  from  the  classical  solutions.  Here,  while  holding  the  mass 
distribution  and  the  moments  of  inertia  of  the  molecule  fixed,  a  branch  of  non-great-circle 
relative  equilibria  starting  with  |A|  =  500  and  ending  with  |A|  =  40000  is  obtained.  In 
all  examples,  the  Newton- Raphson  method  is  used  to  find  solutions  in  each  step  of  the 
continuation  method.  With  the  error  criteria  (6.8),  we  revise  the  standard  algorithm  as 
follows. 

ALGORITHM  7.1. 

1.  Give  xo-  Set  k  0.  Specify  e  . 

2.  Find  8x  such  that 

F'(xk)  8x  =  F(x jfc). 

3.  Let  7  =  max{&r,-,  i  =  1,7}.  Evaluate  a,  a.  If  7  <  e  and  0*77  <  1/2,  EXIT. 

4.  Xfc+i  =  Xk  +  8x.  k  =  k  +  1 .  Goto  2. 

I! 
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In  Step  2,  the  singular  value  decomposition  of  F'(xk)  is  found  first.  This  not  only 
takes  care  of  the  problem  of  ill-conditioning,  but  also  provides  us  with  the  smallest  singular 
value  which  leads  to  the  constant  a ,  cf.  (6.10).  The  other  constant  a  can  also  be  obtained 
from  (6.10).  As  guaranteed  by  Theorem  6.1,  if  Algorithm  7.1  successfully  converges  to 
an  x ,  there  must  exist  a  solution  x *  nearby.  The  error  bound  is  given  by  e  ,  which  will 
be  chosen  to  be  10~8 .  In  the  following  examples,  every  solution  is  obtained  through  this 
Algorithm,  and  is  thus  a  bona-fide  approximate  solution. 


7.1.  Phobos 

We  consider  the  case  of  the  satellite  Phobos  moving  around  the  planet  Mars. 
The  related  data  for  our  study  is  included  in  Appendix  A,  cf.  [15].  In  the  process  of 
nondimensionalization,  the  scales  are  chosen  to  be,  cf.  (2.14), 


mass  :  m  =  1.082  x  1016  (kg),  length  : 


12.4  (km). 


(7.1) 


The  time  scale  is  irrelevant  here,  since  we  are  interested  in  relative  equilibria.  Accordingly, 
the  scaled  radius  is  about  |A|  =  760.  Also  the  non-dimensional  moments  of  inertia  are 


/  0.3294  0  0  \ 

I  =  0  0.2825  0 

\  0  0  0.3881  J 


(7.2) 


Now  we  assume  that  Phobos  is  in  the  shape  of  an  ellipsoid  with  the  above  inertia  and 
uniform  mass  distribution.  The  axes  of  the  ellipsoid  can  be  found  through  elementary 
analysis  to  be  1.847,  2.086,  1.496.  We  next  look  for  the  most  asymmetric  molecule  fitting 
inside  the  ellipsoid  (cf.  Figure  7.1)  to  get  the  most  interesting  non-great-circle  relative 
equilibria.  Physically  this  procedure  might  be  interpreted  as  seeking  the  largest  possible 
orbit  deviation  due  to  non-uniform  mass  distribution  consistent  with  an  observed  shape. 
The  center  of  molecule  is  set  to  coincide  with  the  center  of  the  ellipsoid.  The  constants  a , 
a  for  this  example  are  approximately  104 ,  and  103  respectively.  It  was  found  that  with 

x  1  =  0.9236,  x2  =  0.4618,  ar3  =  1.043, 

x4  =  0.5214,  x5  =  0.748,  x&  =  0.748, 
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Figure  r/  l.  Molecule  inside  the  ellipsoid  of  Phobos 

we  obtain  the  following  non-great-circle  relative  equilibrium, 

A  =  (759.999,1.215,-1.3  x  10“10), 

a  =  (2.6  x  10~n, 0.000, 152.0),  /3  =  -577600.40. 

In  spherical  coordinates  ( r,6,<p ),  the  vector  A  =  (rcos<£cos0,  rcos<£sin0,  rsin^)  can  be 
represented  as  (760.00,0.0916,0.000).  Here  9  and  <p  are  in  degrees.  Clearly  this  relative 
equilibrium  is  very  close  to  the  great-circle  case  in  which  A  is  (760,0,0).  In  order  to 
construct  interesting  non-great-circle  relative  equilibria,  we  hold  the  size  of  the  molecule 
and  the  third  moment  of  inertia  (0.3881)  fixed,  while  varying  the  other  two  moments  of 
inertia.  We  obtain  a  family  of  relative  equilibria  as  listed  in  Table  7.2,  where  the  norms 
of  A  and  fi  are  760  and  152,  respectively.  The  notations  9(\)  and  <j>( A)  denote  spherical 
coordinates  of  the  vector  A  (in  degrees). 

Figure  7.3  shows  the  configuration  of  the  molecule  at  relative  equilibrium  for  the  case 
of  I\  =  0.305956,  I2  =  0.305935  in  Table  7.2.  The  molecule  rotates  approximately  40.13 
degrees  about  the  vertical  axis.  Here  we  see  the  effects  of  nearly  equal  moments  of  inertia 
in  changing  the  configurations  at  relative  equilibria.  It  is  observed  in  Table  7.2  that  as 
the  moments  of  inertia  approach  each  other,  the  molecule  at  relative  equilibrium  exhibits 
significant  change  in  orientation  from  the  classical  solutions. 

We  remark  that,  in  [29],  it  was  shown  that  there  is  no  non-great-circle  relative 
equilibrium  for  the  second-order  approximate  model  with  large  |A| .  Here  the  second-order 
approximate  model  refers  to  the  model  with  the  approximation  of  the  potential  function 
V(A)  defined  in  (3.10)  as 
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h 

h 

0(A) 

0(0) 

m 

0.329386 

0.282505 

0.0916 

0.0000 

0.0509 

90.000 

0.325386 

0.286505 

0.1094 

0.0000 

0.0675 

90.000 

0.320086 

0.291805 

0.1485 

0.0000 

0.1049 

90.000 

0.315186 

0.296705 

0.2251 

0.0000 

0.1796 

90.000 

0.310036 

0.301855 

0.5087 

0.0000 

0.4604 

90.000 

0.308036 

0.303855 

1.0152 

0.0000 

0.9648 

90.000 

0.306886 

0.305005 

2.3956 

0.0000 

2.3414 

90.000 

0.306636 

0.305255 

3.3944 

0.0000 

3.3379 

90.000 

0.306436 

0.305455 

5.0607 

0.0000 

5.0010 

90.000 

0.306336 

0.305555 

6.6530 

0.0000 

6.5906 

90.000 

0.306236 

0.305655 

9.5089 

0.0000 

9.4431 

90.000 

0.306131 

0.305760 

15.7147 

0.0000 

15.6474 

90.000 

0.306086 

0.305805 

20.2084 

0.0000 

20.1450 

90.000 

0.306066 

0.305825 

22.6487 

0.0000 

22.5890 

90.000 

0.306016 

0.305875 

29.8738 

0.0000 

29.8313 

90.000 

0.305996 

0.305895 

33.1418 

0.0000 

33.1095 

90.000 

0.305956 

0.305935 

40.1303 

0.0000 

40.1231 

90.000 

Table  7.2.  Relative  Equilibria  for  Molecules  of  Different  Inertias  at  Constant  Orbit  Radius  |A| 

~  ~\X\  ~  2|Ap  +  2|Ap  <  A’  IA  >  ' 

In  order  to  have  this  second-order  model  fail  to  be  a  valid  approximation,  we  expect  to  have 
inertias  close  to  each  other.  Of  course,  inertial  spherical  symmetry  does  not  necessarily 
imply  a  spherically  symmetric  body,  so  that  the  body  may  be  subjected  to  higher-order 
terms,  cf.  [14].  Our  numerical  investigation  is  in  accord  with  the  above  observcitions. 

In  [1],  small  deviations  from  great-circle  relative  equilibria  were  found  as  perturba¬ 
tions  of  solutions  to  the  approximate  model.  The  numerical  results  presented  here  show 
that  with  two  nearly  equal  moments  of  inertia,  the  exact  problem  may  have  non-great- 
circle  relative  equilibria  with  large  changes  in  orientation  from  the  solutions  of  the  classic 
approximate  model. 


V(A) 
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Figure  7.3.  Configuration  of  a  Non- great- circle  Relative  Equilibrium 


7.2.  Example  with  Branches  of  Solutions  Parametrized  by  |A| 

We  consider  a  molecule  with  the  following  data, 

m1  =  0.330066,  m2  =  0.00330033, 
m3  =  0.330033,  m4  =  0.00330033, 

(7.3) 

m5  =  0.33,  7ng  =  0.00330033, 
h  =  0.3332,  I2  =  0.3335,  I3  =  0.3333, 


moving  in  a  central  gravitational  field.  The  intuitive  idea  behind  devising  example  (7.3)  is 
to  choose  the  moments  of  inertia  close  to  each  other  so  as  to  emphasize  the  importance  of 
higher-order  terms,  while  the  mass  distribution  is  designed  such  that  the  body  is  structurely 
highly  asymmetric.  In  particular,  the  ratios  of  masses  are  101, 100,  and  99  for  each  principal 
axis,  respectively. 

The  method  of  continuation  in  mass  discussed  in  Section  5  has  been  used  with 
Algorithm  7.1  to  find  relative  equilibria  for  selected  |A|  from  500  to  40,000.  We  remark 
here  that  the  nondimensional  parameter  regime  overlaps  the  range  in  which  Earth  orbiting 
artificial  satellites  may  lie.  For  example,  the  space  station  Freedom  has  been  designed  to 
have  maximum  dimension  0.2  (km)  and  orbit  with  radius  6,800  (km).  Thus  the  scaled  |A|  is 
approximately  34,000,  which  falls  within  the  examined  range.  With  some  experimentation, 
it  was  observed  that  with  the  initial  point, 
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A  =  ^|A|j  ,  a  =  (V)  ,  fi  =  -|A|2,  (7.4) 

in  the  method  of  continuation  in  mass,  the  most  interesting  solutions  were  obtained. 
However,  in  finding  solutions  for  large  |A|,  for  example  |A|  =  20,000,  the  constants  are 
approximately 


a  ~  106,  a  ~  105, 

which  implies  that  we  need  to  have  the  constant  7  <  10~12  to  fulfill  the  condition 
acr 7  <  0.5.  Because  for  these  solutions  |A|  ~  105,  the  16-digit  accuracy  available  in 
doable  precision  on  a  typical  32-bit  machine  (e.g.  a  Sun  Workstation)  is  inadequate.  We 
therefore  ran  our  program  on  a  CRAY  supercomputer  which  provides  29  significant  digits 
in  double  precision,  cf.  e.g.  [20].  Table  7.4  lists  those  relative  equilibria  thereby  obtained 
with  A ,  Q,  represented  in  their  spherical  coordinates. 


|A| 

0(A) 

l«l 

6(0) 

m) 

500 

46.8611 

-17.4627 

100 

-54.7456 

-32.6009 

760 

47.8276 

-17.8208 

152 

-54.7761 

-34.1683 

1000 

48.7091 

-18.1277 

200 

-54.8384 

-35.5845 

3000 

55.3232 

-17.2751 

600 

-20.2429 

38.7129 

6000 

65.6627 

-19.5986 

1200 

-15.6572 

22.9702 

8000 

71.0045 

-22.9513 

1600 

-11.4520 

17.2236 

9000 

73.0742 

-25.3278 

1800 

-9.5049 

15.2639 

12000 

78.0382 

22.8848 

.  2400 

-7.5808 

-10.2578 

15000 

81.2009 

18.9175 

3000 

-6.4284 

-6.8820 

20000 

84.1137 

13.7814 

4000 

-4.9693 

-3.7331 

30000 

86.5362 

8.1721 

6000 

-3.2604 

-1.4159 

40000 

87.5514 

5.6183 

8000 

-2.3792 

-0.7051 

Table  7-4-  Relative  Equilibria  Illustrating  Large  Orientation  Drifts 
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Interestingly,  there  are  three  segments  of  distinct  branches  of  solutions  shown  in 
Table  7.4,  namely,  those  with  |A|  in  the  range  of  500-1,000,  3,000-9,000,  and  12,000- 
40,000.  We  conjecture  that  there  axe  some  bifurcation  or  turning  points  in  the  (A,  j3) 

space.  In  order  to  construct  one  branch  completely,  we  abandoned  continuation  in  the 
masses,  and  instead  started  with  solutions  at  |A|  =  12,000  and  reduced  the  radius  | A J 
with  numerical  continuation  until  ] A]  =  500.  On  the  CRAY,  this  method  of  continuation 
in  radius  successfully  gives  us  the  branch  of  solutions  listed  in  Table  7.5.  It  can  be  seen 
that  the  solution  branch  for  |A|  =  500-1000  in  Table  7.4  is  disconnected  from  the  branch 
starting  from  |A|  =  12,000-40,000  that  is  continued  in  Table  7.5  down  to  small  radii. 


I'M 

0(A) 

<K  A) 

m 

e(a) 

m 

500 

46.7440 

35.1556 

100 

-10.1838 

-37.7702 

1000 

48.5200 

35.0055 

200 

-10.3029 

-36.4710 

2000 

52.0762 

34.5767 

400 

-10.4690 

-33.7789 

3000 

55.5958 

33.9734 

600 

-10.5296 

-30.9911 

4000 

59.0231 

33.1958 

800 

-10.4787 

-28.1563 

5000 

62.3003 

32.2504 

1000 

-10.3185 

-25.3351 

6000 

65.3718 

31.1514 

1200 

-10.0596 

-22.5938 

7000 

68.1912 

29.9214 

1400 

-9.7205 

-19.9952 

8000 

70.7275 

28.5895 

1600 

-9.3239 

-17.5892 

9000 

72.9684 

27.1885 

1800 

-8.8934 

-15.4074 

10000 

74.9195 

25.7512 

2000 

-8.4496 

-13.4624 

11000 

76.6003 

24.3081 

2200 

-8.0085 

-11.7507 

12000 

78.0382 

22.8848 

2400 

-7.5808 

-10.2578 

15000 

81.2009 

18.9175 

3000 

-6.4284 

-6.8820 

20000 

84.1137 

13.7814 

4000 

-4.9693 

-3.7331 

25000 

85.6324 

5000 

-3.9630 

-2.2060 

30000 

86.5362 

8.1721 

6000 

-3.2604 

-1.4159 

34000 

87.0289 

6.9333 

6800 

-2.8441 

-1.0438 

35000 

87.1310 

6.6760 

7000 

-2.7551 

-0.9729 

40000 

87.5514 

5.6183 

8000 

-2.3792 

-0.7051 

Table  7.5.  A  Branch  of  Non- great- circle  Relative  Equilibria 
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It  would  be  interesting  to  find  the  other  branches  completely.  However,  more  delicate 
numerical  tools  need  to  be  applied.  In  particular,  a  continuation  code  designed  to  traverse 
turning  points  is  apparently  required.  This  problem  is  currently  under  investigation. 


8.  Conclusions 

In  this  paper,  we  have  proven  the  existence  of  non-great-circle  relative  equilibria  for 
the  exact  model  of  motion  of  a  rigid  body  in  a  central  gravitational  field.  A  coordinate- 
free  approach  was  adopted  to  avoid  singularities  and  cumbersome  manipulations  inherent 
to  local  coordinates,  such  as  Euler  angles.  Numerical  methods  for  finding  the  relative 
equilibria  with  an  associated  delicate  error  analysis  have  been  presented.  The  numerical 
solutions  we  obtained  demonstrate  that  the  non-great-circle  relative  equilibria  may  exhibit 
large  orientation  changes  from  the  classic  approximate  solutions  for  bodies  with  nearly 
equal  moments  of  inertia.  If  these  relative  equilibria  prove  to  be  stable,  the  attitude  design 
of  spacecraft  might  be  able  to  take  advantage  of  such  configurations  to  yield  effective 
attitude  control  strategies. 
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Appendix  A.  Data  for  Phobos  around  Mars 

Phobos  mass  ~  1.082  x  1016  (kg). 

Phobos  volume  ~  5.673  x  103  (km3). 

Phobos  principal  moments  of  inertia 

/  5.50  0  0  \ 

r  ~  0  4.718  0  x  1017  (kg-  km2). 

\  0  0  6.481  / 

Phobos  orbit  radius  (major  semi-axis)  ~  9,378.5  (km). 

Phobos  radius  ~  11.1  (km). 


30 


Mars  radius  ~  3,394  (km). 


Appendix  B.  Determining  Constants  in  Theorem  5.1. 

In  this  Appendix,  we  discuss  the  issue  of  determining  the  constants  in  Theorem  6.1 
as  well  as  how  to  pick  a  suitable  constant  A  to  make  the  problem  numerically  tractable. 
We  first  look  at  cr,  the  constant  associated  with  the  Jacobian  of  F.  Since  (6.1)  comes 
from  a  variational  principle,  the  derivative  of  F  in  (6.5)  with  respect  to  A,  Q,  and  (3  is  a 
7x7  symmetric  matrix  with  real  eigenvalues.  Let  pm\a  denote  the  minimum  eigenvalue 
of  F'(A,f2,/3).  It  is  easy  to  see  that  l/pmin  is  the  maximum  eigenvalue  of  F'(A,  Cl,  /3)-1  • 
Thus 

||f'(A,n,J8)-'|U  <  — ,  (B.l) 

P  min 

where  ||  •  ||2  denotes  the  matrix  or  operator  norm  induced  from  the  usual  vector  2-norm 
on  IR" .  Accordingly,  if  the  2-norm  is  used,  a  could  be  specified  as  1/ pm\n  •  On  the  other 
hand,  if  the  sup-norm  defined  by  Uxll,*,  =  max{|x,-|,i  =  l,n)  is  selected,  because  of  the 
inequalities 

IMloo  <  ||x||2  <  nllxHoo,  for  x  6  IRn, 


we  have 

lir(A,n,«-‘iu  <  nlir'tA,  $!,«->  ib  <  — . 

P  min 

The  constant  a  can  then  be  chosen  as 

n 

a  =  - . 

Pmin 


(B.2) 


(B- 3) 


Since  it  is  considerably  more  convenient  to  use  the  sup-norm  in  determining  the  Lipschitz 
constant  for  F' ,  the  a  defined  in  (B.3)  will  be  adopted. 

In  order  that  conditions  (6.8)  can  hold,  it  is  evident  that  a  cannot  be  too  large,  or 
equivalently,  pmin  cannot  be  too  small.  In  order  to  have  a  measure  of  pmin ,  we  consider 
the  special  case  of  the  relative  equilibrium  (6.2)  for  a  symmetric  molecule. 
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The  Jacobian  of  F  can  be  found  from  (6.1)  to  be  the  symmetric  matrix, 


+  EUj^fv^  +  Qi)(^Qi)T 

V 


-(A  •  0)1  +  2A0T  -  0AT 
I  +  (I'M2  +  ^)1  -  AAT 


°\ 

0 

0/ 

{BA) 


Evaluation  of  F'  at  the  particular  relative  equilibrium  (6.2)  and  (6.4)  yields  the  following 
simplified  matrix, 


/  Jn 

0 

0 

0 

2  ab 

0 

°\ 

0 

J22 

0 

—ab 

0 

0 

0 

0 

0 

Jz3 

0 

0 

0 

0 

0 

—ab 

0 

h-h-  a 2 

0 

0 

0 

(5.5) 

2  ab 

0 

0 

0 

0 

0 

b 

0 

0 

0 

0 

0 

h~h 

0  j 

V  0 

0 

0 

0 

b 

0 

0  / 

where 


Jn  =  b2  +  A 


3  mi 


+ 


3mi 


+ 


6m3<i2 


+ 


6m5a2 


mi 


(a  +  Xly  (a-xO3  ^a2  +  xf 


J, 


6m3X3 


22 


-E 


\J  a2  +  xf  ^\X  +  Qi\3J' 


J33  =  b2  +  A 


6  ra5x| 
y/a2  +  x|5  + 


-E 


IriA+Qii3;’ 

(5.6<z) 

(B.6b) 

(5.6c) 


The  chaxacteristic  polynomial  of  J  is 

det(sl-J)  =  (5  -  J33)(s  -  (I3  -  J2)) 

(s2  —  (J22  +  h  —  h  —  a2)s  +  J22i.I1  —  I2  —  a2)  —  a2b~) 
(s3  —  Jus 2  —  (1  +  4  a2)b2s  +  Jn&2^ . 


There  are  two  obvious  eigenvalues  of  J ,  namely,  J33  and  h  —  h-  The  second  is  given  by 
body  parameters,  and,  in  our  nondimensionalization,  has  absolute  value  smaller  them  one. 
With  (6.7),  J33  in  (B.6c)  can  be  approximated  as 
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A  6Am5xl 
J33  ~  -=•  + 


E 


Ami 


y/a2  +  xf  ,-=ilA  +  (5>l 


§Am*,x\ 

y/a?  +  xf 


+4|>G?-^W) 


Aim  1 


xx[(g  +  xi)2  +  g(g  +  Ji)  +  a2]  m ^  xi[(g  -  xi)2  +  a(o  -  x\)  +  a‘1}  (B.8) 


g3(g  +  xi)3 


g3(g  —  X\ )3 


+  2m3 


(3/g2  +  x 3  —  g)[g2  +  rcg  +  g>/g2  +  a:2  +  g2 

g3  ^/g2  +  X3 


( y/ g2  +  xf  -  g)[g2  +  x^  -f  g-^/a2  +  xf  +  g2]  6m5i5 

*4“  Zms  5  1 


g3  4/g2  +  x§ 


V«^  + 


x; 


It  is  now  clear  that  the  magnitude  of  J33  is  at  most  A/g4.  Thus  in  order  to  have  J33 
large,  a  reasonable  criterion  for  A  is 

A  ~  g4.  (5.9) 


Next  we  check  the  third-order  term  in  (B.7).  The  coefficient  for  s  is  approximately 

(l+4a2)62  ~  (l  +  4a2)4  =*  — ,  (B- 10) 

v  a3  a 

for  large  a.  The  choice  of  A  in  (B.9)  makes  this  coefficient  large  which  in  turn  implies 
there  is  a  large  eigenvalue.  As  a  consequence,  the  condition  number  of  the  J acobian  matrix 
becomes  large  for  large  g,  or  large  radius  (with  the  size  of  molecule  fixed).  This  illustrates 
that  the  problem  of  solving  (6.1)  is  highly  ill-conditioned.  However,  by  using  Theorem  6.1, 
this  difficulty  can  be  handled. 

We  next  check  the  other  constants  in  Theorem  6.1.  By  finding  Sx  €  IR,"  such  that 

F'(xq)6x  =  F(x  0),  (5.11) 


the  constant  7  can  be  specified  as 

7  =  ||6x||oo  =  max{<§xj,  i  =  1, •••,«}•  (5.12) 

However,  in  solving  (B.ll)  for  Sx,  since  F'(x 0)  is  ill-conditioned,  we  need  to  use  a  robust 
algorithm.  A  method  based  on  Singular- Value  Decomposition  (SVD)  is  chosen  to  deal 
with  such  difficulties. 
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The  remaining  constant  a  is  the  Lipschitz  constant  for  F' .  Since  F'  is  continuously 

differentiable  at  A  ^  0,  we  first  find  a  bound  on  \\F"\\ ,  where  F" ,  a  third-order  tensor,  is 

the  derivative  of  F'  with  respect  to  (A,f2,/3).  Moreover,  F"  can  be  represented  as  a  linear 

2 

operator  from  lRn  to  1R"  ,  which,  in  turn,  can  be  expressed  as  an  (n2  x  n )  matrix  P .  It 
is  easy  to  show  that 

IIFloo  <  n max{P,j,  i  =  l,n2,  j  =  l,n}.  (I?. 13) 


Because  of  the  simplicity  of  this  bound,  the  sup-norm  is  adopted  in  Theorem  6.1.  Since  in 
the  Theorem,  the  constant  a  should  be  a  Lipschitz  constant  on  the  whole  domain  D ,  we 
look  for  a  symbolic  representation  of  a .  It  can  be  checked  that  the  elements  in  the  matrix 


P  have  standard  forms:  2A ,  2f2 ,  and 


E—  l5Arm\\  .  .  6Am,-A3  ZArrii 

•  Q>)  +  L  '  .  „  +  L—  at;(a  +  9.) 


15Am5(A3  +  rr5)2 
|A  +  Q5|7 
15Am6(A3  —  x6)2 
1A  + Qel7 


fe|A  +  <fc| 


|A  +  Qi\ 


t  \  ,  6Am5(A3  +0:5) 

(A+<?5)  + . iu'w  63 

,  6Am6(A3  -  x6) 

(X+Qe)  +  <v 


where  e3  =  (0  0  1)T .  Let  mmax  =  max{mj,  i  =  1,  •  •  • ,  6} .  Bounds  for  these  terms  can  be 
found  to  be 


{  2|A|,  2|ft|, 


144Amr 

"IF 


(.B.14) 


In  order  to  monitor  the  bound  conveniently,  we  select  A  such  that 

144Ammax 


\M* 


<  2|A|, 


which  is  equivalent  to 


A  < 


1A|5 

7277lmax 


(B.15) 


This  constraint  on  A  associated  with  condition  (B.9)  provide  us  a  good  criterion  to 
determine  A  such  that  the  numerical  method  is  tractable. 

With  (B.15)  and  (6.7),  we  can  approximate  |Q|  as 
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\n\  ~ 


A 

Tap 


1/2 


(B.1G) 


< 


|A| 


y/n 


mr 


<  |A|, 


since  mmax  >  1/6.  Thus  (B.16)  and  (B.15)  ensure  that  the  maximum  value  in  P  is  2|A|. 
A  symbolic  bound  of  a  can  then  be  chosen  to  be 


a  ~  2n|A|. 


(BA7) 
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